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AN ANALYSIS OF THE TRANSIENT SOLIDIFICATION OF A FLOWING WARM LIQUID ON A CONVECTIVELY COOLED WALL 

by Robert Siegel and Joseph M. Savino 
Lewis Research Center 

National Aeronautics and Space Administration 
Cleveland, Ohio 

ABSTRACT 

A study is made of the frozen layer that forms when a warm liquid flows over a flat plate that is 
cooled below the freezing temperature of the liquid by a coolant flowing along the other side of the plate 
Three analytical procedures are employed and compared for accuracy and convenience in application. One is 
an iterative procedure, yielding three analytical closed-form solutions, each successive solution being a 
more accurate, higher order approximation. The two other analytical procedures require numerical integra- 
tion for evaluation. Numerical results are presented graphically so that rapid estimates of frozen layer 
growth can be made. 


AUSZUG 

Die Studie befasst sich mit der erstarrten Schicht die entsteht wenn eine warme Fliissigkeit uber 
eine flache Platte fliesst und diese Platte durch ein Kuhlmittel an der Unterseite auf den Gefrierpunkt 
der betreffenden Flussigkeit gebracht wird. Die drei verwendeten analytischen Methoden werden in Bezug 
auf Genauigkeit und Bequemlichkeit der Anwendung verglichen. Eine von diesen drei Methoden ist iterativ 
und ^liefert drei analytische Losungen in geschlossener Form, in denen jede Folgende eine genauere 
Annaherung einer hoheren Grossenordnung ist. Die beiden anderen analytischen Methoden erfordern 
numerische Integration zur Auswertung. Numerische Ergebnisse werden graphisch dargestellt und machen da- 
her rasches Abschatzen des Wachstums der erstarrten Schicht moglich. 


AHHOTAUMfl 

UccjieayeTCH 3aMopoaceHHbiH cjioh o6pa3yiomiificH npn oGTeicaHHH iijiocko 8 njiacTHHKH Teiuiog 
sch^koctbio c o,h;ho 8 cTopoHhi h oxjia^HTejieM c npyroft ctopohm. Hcnojib3yioTc,s Tpn aHajiHTHwecKHX 
npneMa KOTopwe cpaBHHBaiOTCH no towhoctk h y^oGcTBy npnjiosceHHH. O^hh HBJiHeTCH HTepaTHBHUM 
npneMOM flaiomuM Tpn 3aMKHyTbix pemeHHft Kasc^oe H3 kotophx npe^cTaBJineT Goiiee TOWHoe npn6jin- 
sceHHe Bhicmero nopn^Ka. OcTaiomHecn £Ba aHajiHTHwecKHx npneMa aojukhu oneHHBaTcn wucjieHHhiM 
HHTerpnpoBaHHeM. ^HCJieHHiJe pe3yjibTaTbi npe^CTaBJieHH rpa$HwecKH GuCTpofi oneHKH pocTa 
3aM0p03R6HH0r0 CJIOH. 


INTRODUCTION 

This paper is concerned with the transient 
solidification of a warm flowing liquid that is in 
contact with a wall that is cooled below the freez- 
ing point of the liquid by a coolant flowing on 
the opposite side. This type of solidification 
problem arises in important and familiar applica- 
tions such as continuous casting of metals and 
freezing of rivers. The application that motivated 
the present study is concerned with a number of 
advanced propulsion devices that utilize liquid- 
liquid heat exchangers in which the coolant is at a 
temperature below the freezing point of the warm 
liquid. If the temperature or flow rate of the 
warm liquid is low enough, a frozen layer will form 
on the exchanger walls. The prevention of serious 
flow blockage by excessive solidification requires 
an understanding of the transient formation of the 
solid layer on the exchanger walls. 

The phenomenon of freezing (or melting) of a 
material has been the subject of considerable math- 
ematical analysis for two reasons. The first is 
because of its practical importance in applications 
such as mentioned above. The second is that the 
problem offers a challenge for mathematical study 
because the liquid-solid interface is a moving 
boundary, which provides a nonlinear mathematical 
condition. Exact analytical solutions for the 
transient frozen layer thicknesses and heat trans- 
fer are usually difficult to obtain except for 
boundary conditions of certain types. This has led 


to the situation where in some instances physically 
unrealistic boundary conditions have been imposed 
in order to achieve a mathematical solution. Some 
exact solutions for conditions of practical impor- 
tance are summarized by Car slaw and Jaeger [l] 1 . 

No attempt will be made here to review the mathe- 
matical bibliography which is not directly perti- 
nent to the present study. A few references which 
provide various mathematical melhods and have lists 
of references are given by [2-5]. 

When an exact analytical solution cannot be 
found, numerical or approximate methods must be em- 
ployed. Murray and Landis [6] devised a numerical 
finite difference technique and demonstrated it for 
two sample problems. Cochran [7] developed a 
lumped -parameter method, wherein all the heat ca- 
pacity is lumped at the center of the frozen mate- 
rial. A heat balance integral approach was demon- 
strated on a variety of freezing and melting prob- 
lems by Goodman [8] where the temperature distribu- 
tion in the material was approximated by a second 
or third order polynomial to evaluate the integrals 
of the conduction equation. Adams [9] employed an 
iterative technique to find the frozen layer thick- 
ness and heat transfer in metal castings with and 
without accounting for the conduction into the 
mold. This method does not involve any approxima- 
tions in principle, but provides an approximate re- 
sult because in some instances only the first 
iteration can be conveniently carried out in 
analytical form. 

^Number in brackets denote references. 
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The model chosen for study here is a one- 
dimensional solidification of a liquid that flows 
over one side of a thin flat wall (Fig. 1). The 
other side of the wall is being convectively cooled 
by a fluid that is at a temperature below the 
freezing point of the liquid. The wall on which 
the liquid solidifies is assumed to be sufficiently 
thin so that the heat flow needed to subcool the 
wall during the transient is negligible compared 
with the heat flow through the wall. This heat flow 
is comprised of the latent heat of fusion, the sub- 
cooling heat capacity of the frozen layer, and the 
convection at the solidified interface. 

Three analytical methods will be employed. For 
two of these the transient heat conduction equation 
which governs the energy balance within the solid- 
ified layer has been integrated in a general way to 
provide an expression for the frozen layer growth as 
a function of time. The general integrated equation 
is evaluated for the present problem by using a 
method of successive analytically iterated approxi- 
mations and also by using an approximate temperature 
profile of the type proposed by Goodman [8]. The 
former technique has been utilized in a different 
manner by Adams [9] for some other problems. The 
analytical approximations are carried out here in an 
improved manner so that higher order solutions can 
be found, and an analytical expression very close to 
an exact solution is obtained. The third analytical 
method is one proposed by Goodman [8] where a heat 
balance at the frozen layer -liquid interface is 
used. The methods are compared with regard to accu- 
racy and ease of application. One of the closed 
form solutions represents for all practical purposes 
the exact solution and is used to evaluate the final 
results that are presented in graphical form. 


ANALYSIS 


The model chosen for this study is the one- 
diraensional configuration shown in Figure 1. A warm 
liquid at a fixed temperature Tj flows over one 
side of a thin plane wall providing a constant con- 
vective heat transfer coefficient hj. A coolant at 
a fixed temperature T c flows over the opposite 
side of the plate providing a constant convective 
heat transfer coefficient h c . A transient solidi- 
fication process can then be initiated in a number 
of ways such as by reducing the temperature or the 
flow rate of the warm liquid, introducing the cool- 
ant if it had not been flowing, or lowering T c of 
the flowing coolant. These changes will cause the 
wall to further cool until the freezing temperature 
is reached on the surface of the plate exposed to 
the warm liquid. At this instant (x = 0, X(x) = 0), 
since the specific heat of the wall has been ne- 
glected, a linear temperature distribution will 
exist in the plate as shown in Figure 1, and solid- 
ification is assumed to begin. 

As the solid layer on the wall increases in 
thickness, heat is being extracted from the boundary 
layer of the warm liquid and is transferred by con- 
vection to the liquid-solid interface. There it 
combines with the latent heat released by the solid- 
ification process and is conducted through both the 
frozen layer and cold wall to the wall surface in 
contact with the coolant, where it is convectively 
removed. An additional amount of heat removal is 
needed to subcool the frozen layer, and this is also 
transferred to the coolant (the heat removal neces- 
sary to subcool the wall is neglected in comparison 
with the other heat flows). The solid layer con- 
tinues to grow until it achieves a steady- state 


thickness X s . Constant properties are assumed 
throughout the analysis. 


Steady- State Thickness of Frozen Layer 

In many solidification problems the frozen 
layer never approaches a steady- state thickness ; 
rather it continues to propagate with time into the 
liquid. When, however, the liquid flowing over the 
frozen layer is at a temperature above the freezing 
point, the frozen layer will achieve a steady-state 
thickness. This thickness will be used as a refer- 
ence length in the later analyses of growing layers 
and hence will be separately derived here. 

It is assumed that the liquid-solid interface 
at the boundary of the frozen layer is at the 
freezing temperature Tf. If the heat flow is 
taken as positive in the positive x direction, 
the convected flux from the liquid to the liquid - 
solid interface is 


-q = hj(Tj - T f ) = constant 


Using the overall heat transfer resistance between 
the liquid-solid interface and the coolant gives 
the relation 


-q = hj(T , - T f ) 


This is solved for the steady-state thickness 


^ = 


k T f - T c 

hi T : - Tf 



( 1 ) 


For X s = 0, Equation (l) gives the relation be- 
tween variables required to just avoid freezing. 
For example, solving for T i gives 



In order to prevent freezing the liquid temperature 
Tj must be equal to or greater than the value 
given by Equation (2). 


General Equations for Frozen Tayer Growth and 
Temperature Distribution 


The heat flow within the frozen layer is 
governed by the transient heat conduction equation, 


. S 2 T St 

k y> = pcp ^ 


(3) 


Equation (3) is integrated from any position x 
within the layer to the solid-liquid interface X 


ST 

37 




St 

37 


dx 


(4) 


At the interface the heat conducted into the solid- 
ified layer is equal to that supplied by the latent 
heat of fusion and the convection from the flowing 
liquid. 


k s| = pL if + h * (T i ' Tf) 


(5) 
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Equation (5) is substituted into Equation (4) to 
give 

' k !| x = -^37- h * (T * - ^ + pc p/ X I 

(6) 

The term on the left side is the heat flow crossing 
any position x at any time t. The last term on 
the right is the heat removed to subcool the portion 
of the solidified layer between x and X. The in- 
tegration of Equation (6) from the wall (x = 0) to 
any position x results in an expression for the 
instantaneous temperature distribution 

T(x,r) - T 2 = £ H x + h. (Tj - T f )x 


£f] 

k 


l [ L * “H 


( 7 ) 


_ 1 1 Rx' /dX' \ RS [l / X ST' 
■ 1 + R ^dT' +1 y + i+RR < yr ) ST 7- 


dx' 


cu: 


( 10 ) 


As a consequence of the nondimens ionalizat ion it is 
noted that two parameters S and R have appeared. 
Their physical significance is discussed in the 
section RESULTS AND DISCUSSION. 

By applying the rules for differentiating 
under an integral , the integral terms of Equa- 
tion (10) can be transformed so that the derivative 
is taken outside the integral signs, and Equa- 
tion (10) becomes 


T’ = 1 - 


1 + Rx' /ax' \ RS b (l r X '. 

^mc\6^ + V - ITr 57^ r 


dx' 


tion of time, Equation (6) is written at x = 0 to 
give the heat flow at the wall, and this is equated 
to the heat flow through the wall and into the cool- 
ant, 


q(x = 0) = -pL S - h,(T, - T f ) 

dr ‘ ‘ 1 


P C P 


/ 


3T aT _ T 2 ‘ T C 

S 7 1 7 a 


(a) 


r/j 


T' dx' dx 


( 11 ) 


By reversing the order of integration, the double 
integral in Equation (11 ) can be transformed into 
equivalent single integrals. This is accomplished 
by first writing the double integral in terms of 
dummy variables of integration q and p as 
follows : 


Then Equation (8) is solved for Tg and substituted 
in Equation (7) to give 




pL 


+ hj(T z -T f )-pc 


l H* 






dxjdx (9) 


Equation (9) is then placed in dimensionless form 

■by letting x' = x/X s , X' - X/X s , T ’ = hj(T z - T f )r/ ^ x . - X ' 

pLX s , and T' = (T - q)/(T 0 - Tf) and using Equa- /, T (x ,X ) dx dx' 

tion (1) in the form (X g h^/k) (Tj - T f )/ 

(Tp - T . ) = R/(l + R) to cive 


f x' ,X’ 

/ / T 1 (3,x‘ ) dp dn 


where X* is treated as a constant. Reversing the 
order of integration leads to 


r x r x r x ' r 0 r*' r *' 

/ / T'(P,X') dp dn = / / T'(p,X') dq dp + / / 

J 0 * / n y O J D Jr' Jn 


T'(p,X') dq dp 


/ x 1 /*X ’ 

pT’(P,X’) dp + / x’T'(p,X' ) dp 

x' 


so that 


f x r x r x ' r x ' 

J J T*(x' ,X’ ) dx 1 dx 1 — / xTtx'jX') dx f + x’ / T'fxSX') dx' 

0 •V *'0 J*' 


( 12 ) 


Equation (ll) then takes the form 


T ' (x ' ,X' ) = 1 - — f T' dx* + x 1 f T’ dx* + f x’T 1 dx’) (13) 

1+R \dT / 1 + R Sx' \R Jq x' J 0 / 
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For 


convenience let I(x',X') - ^ / T’ dx' 

R J 0 

/*’ f x’ 

+ x’ J T* dx' + x’T' dx' and note that 

8I/8t' = (8l/8X' ) (dX'/dt ' ). As a result, Equa- 
tion (13) takes the form 

r = 1 - ±^3£-(2L +i (i4) 

l+R ^dT ’ y l+R^TdT 1 

When Equation (14) is evaluated at the solid- liquid 
interface x’ = X’, there results 

_ 1 + RX' /dX' \ RS dG(X’) dX’ 

1 + R ^dr' + y + 1 + R bx' dr 7 ^ 15} 

where 

G(X’)s I(x' = X’,X’) = i /' X T'dx' + f* x’T’ dx' 
R *v\ 


rV /.X 

R *'0 

Solving Equation (15) for dX'/di:' gives 
R(1 - X 1 ) 


dX’ 

dx' 


1 + EX' + RS 


dG 

dX' 


( 16 ) 


Separating the variables and integrating yields 

T ’=/ X ' l|T^“ ,+ ^ X, TTfr (17) 

When the first integral on the right is evaluated 
and the second integral is integrated by parts, the 
final expression for t’ = t'(X') is 


(18) 


This is the expression that will be used to compute 
the frozen layer thickness as a function of time. 
Since G(X') contains the temperature distribution, 
this distribution is obtained from Equation (14) 
with Equation (16) used to eliminate dX' /dx ' . The 
final relation for the temperature distribution in 
the frozen layer becomes 



. , /l+Rx’+RS 

T t = RQx 1 ) ( 8x 

l+R \ l+R 



(19) 


Equation (19) is an integral equation for the 
temperature T'(x',X'). To solve Equation (18), 
Equation (19) tnust first be solved so that the 
integrals that* comprise G in Equation (18) can 
be evaluated. Two methods for carrying this out 
will be examinee! here. 

The first method, and ultimately the more ac- 
curate of the two, is an iterative technique where- 
by the integrals in I and G are evaluated using 
temperature distributions that are obtained from 
successively better approximations of Equa- 
tions (19). It appears that this procedure will 
ultimately lead to an exact solution. In the sec- 
ond method Equation (19) is not used. Instead 
the temperature distribution to be used in G(X) of 
Equation (18) is approximated by a second order 
polynomial with coefficients that are evaluated 
using the physical conditions at the boundaries of 
the frozen layer. Since this approximate temper- 
ature profile does not arise from the energy equa- 
tion, the second method cannot yield an exact solu- 
tion. The two methods are outlined in the next 
sections. 

Solution by Analytical Iterations 2 

An iterative procedure has been utilized by 
Adams (9] for sane solidification problems, and 
before discussing the present method, a few com- 
ments are in order. When Adams' procedure was ap- 
plied to the present problem, the first iteration 
(second approximation) resulted in a cumbersome 
quadratic equation for dX'/d-r'. After solving the 
quadratic equation, a numerical integration was 
required to determine X'(t'). A second iteration 
would have been extremely difficult. Hence an al- 
ternate approach was devised in which higher order 
iterations could be analytically obtained thereby 
leading to more accurate results. These will now 
be presented. 

A first order approximation to the temperature 
distribution and growth times is found by neglect- 
ing the effect of heat capacity within the frozen 
layer. When the Cp becomes zero, the parameter 
S = 0, and Equations (18) and (19) reduce to 


l + R 


In (1 - X’) 


, = R(X* - x') 

I 1 + RX* 1 + RX' 


RX 1 


(1 - i) 


( 20 ) 


( 21 ) 


A second approximation for T* and 
t 1 = t'(X') was obtained by substituting Equa- 
tion (21) for Tj in Equations (18) and (19) with 
the S terms now retained. After the integrations 
in I(x',X') and G(X' ) are performed, the results 
can be simplified to 


T ii - T i - s 


f BX ' |~r + - - - 1 + 5 * 58 + r2 in (1 - X*) + ^ (RX ' + 
|3(RX' + 1) [_ R(R + !)J 3(1 + R) 2 3R(1 + R) 2 J 


and 

T . R(i - x’S) 

1 + RX'| + RS 

81 

dX* 

I 

R(1 - X') 


T n " l + R 

l + R 


1 + RX’ + RS §r| 


( 22 ) 


(23) 


This will be referred to later as method "Al" 
(Analytical Iteration). 
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where 



When the procedure is repeated, that is, when Equation (23) is substituted into Equations (18) and (19), 
the third approximation results, which is the final approximation evaluated here, that is, 
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+ RX*) 3 






X I 2R(l - X')(l + RX’) - (1 + RX') - R(1 - X’)(l + RX 1 ) 3R(l + RX' ) 




( 31 ) 


No attempt was made to analytically evaluate the 
integral in Equation (26) as there did not seem to 
he any advantage to doing so because of its alge- 
braic complexity. Since the entire equation was to 
be evaluated on a digital computer, the integration 
was performed numerically. 

It is obvious from the complexity of the tem- 
perature distribution Equations (27) to (31) that a 
fourth analytically determined approximation would 
be almost impossible to carry out. As will be 
shown, a fourth approximate solution is unnecessary 
because the AI method converged rapidly with the 
first three approximate solutions. 


Approximate Solution by Prescribed Temperature 
Method 3 


The basic expression (Eq. (18)) that predicts 
the dimensionless growth times t 1 for a layer of 
thickness X 1 contains integrals, in the function 
G, of the instantaneous temperature distributions 
in the frozen layer. Because the temperature dis- 
tributions appear as integrands, the growth time 
t* is not very sensitive to the exact shape of the 
distribution. Therefore, lasing a reasonable ap- 
proximation for the temperature profile should re- 
sult in fairly accurate predictions for 
t* = t'(X'). This idea was first applied to prob- 
lems of freezing and melting by Goodman [8], 

For the present problem, the temperature pro- 
file within the frozen layer is approximated by a 
second order polynomial of the form 

T' = A(X* - x') + B(X* - x') 2 (32) 

Equation (32) already satisfies the boundary condi- 
tion that the temperature is equal to Tf at the 
liquid-solid interface. Two additional boundary 
conditions are used to evaluate A and B. The 
boundary condition at the wall is that the follow- 
ing heat balance be satisfied: 



which has the dimensionless form 

- R[T'(x’ = 0) - 1] (33) 

x'=0 

Substituting Equation (32) into Equation (33) pro- 
vides the condition 

-A - 2X'B = R(AX‘ + BX' 2 - 1) (34) 

The second boundary condition is at the liquid- 
solid interface and is given by Equation (5). 


dT* 


2 

This will be referred to later as method "PT" 
(Prescribed Temperature). 


Equation (5) is modified in order to express it 
entirely in terms of T, without dX/dr appear- 
ing, by using the following relation which can be 
written for any location in the frozen layer: 


dT 



dx + 



dr 


At the instant when the interface moves through 
position x, x = X, dx/dr becomes dX/dx, and the 
derivative dT/dx = 0 since the interface is al- 
ways at Tf. Then 


r dT 

1 

Sx St 

J dx 


V — 

31 

II 

i 

31 


at interface 



at given timej 



at fixed x as 
interface moves 
across that x 


This is used to eliminate dX/dT from Equation (5) 


dT 

dx 


St 

dx 

m .& — 

k St 

x £ 


h 7 

+ 1 T (T I - T f> 


If the relation dT/dx = (k/pc^) (d 2 T/dx 2 ) from the 

conduction equation is used, this boundary condi- 
tion can be rearranged into 


l St 

\ _ _ L_ d% 

+ — (T- - T ) — 

1 Sx 

x/ " C P 3x2 

k K l V dx 

X 


In dimensionless form it becomes 


Y _ i s 2 t’ 


R dT' 

x-/ S 3x ' 2 

X' 

1 + R dx' 


(35) 


X' 


The temperature profile (Eq. (32)) is substituted 
into Equation (35), and the result is rearranged 
into 

B = I( aZ - FTe) < 36 > 

which gives B when A is known. Equations (34) 
and (36) are solved simultaneously to yield A 

A = r i + rx' + [ r R 

2 ( 1 + R) S(2X* + RX' 2 ) \L 2 ( 1 + R) 

+ f + 2 

2 I 2 

S(2X* + RX’ )J S(2X* + RX' ) 
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The coefficient A is the root of a quadratic 
equation, and the positive sign in front of the 
square root was chosen so that A would approach 
the proper limit as steady state is achieved. The 
limit is found by realizing that the steady-state 
temperature profile in the frozen layer must be 
linear and hence from Equation (32), B(X* -*■ 1) = 0. 
Then frcm Equation (36), A(X* -► 1) = R/(l + R), 
which is also the limit of Equation (37). 

Equations (36) and (37) give A and B as 
explicit functions of X*; hence, the temperature 
distribution. Equation (32), is given as a func- 
tion of X 1 and x*. The temperature distribu- 
tion is then used to evaluate the integrals in 
Equation (18) 

The quantity G(X') is 

while the quantity 

f X ' o(x„',)_ ax- 

J 0 (l-x') 

has to be integrated numerically because the 
and B contained in G are very complicated func- 
tions of X*. Equation (18) is then evaluated to 
give the time required to form a thickness of 
frozen layer, t 1 = t^X* ), for any values of the 
parameters R and S. 

With A and B known as a function of X* 
and hence as a function of t*, the temperature 
profile in the frozen layer can be evaluated at 
any time during the layer growth from Equation 
(32), which can also be placed in the form 

T' = AX'(1 - I) + BX ,2 (1 - |) 2 (40) 


(38) 


(39) 

A 


Determination of Frozen Layer Growth by 
Integrating Heat Balance at Liquid- 

Solid Interface 4 


This alternate approach suggested by 
Goodman [8] begins with the interface condition 
that follows Equation (34) 


dX 

dx 


8T 

St 

St 

dx 


k a 2 T 
p c p 3x 2 
oT 


dx 


X 


This has the dimensionless form 


dX 1 

dT* 


d 2 T* 

N , 2 

1 + R dx 1 
RS dj* 
d?" 


X* 


(41) 


The approximate temperature profile. Equa- 
tion (32), is used to evaluate the derivatives 
at the interface to give 


4 This will be referred to later as method 
"IF" (Interface). 


dX 1 1 + R 2B 
dr' ' RS A 

The quantity B is eliminated by using Equa- 
tion (36), and the variables are then separated. 
The result is integrated to yield 



This was integrated numerically for various R 
and S using A as a function of X* from Equa- 
tion (37). As in the previous method, the tempera- 
ture distribution is found from Equation (40) with 
B obtained from Equation (36). 


RESULTS AND DISCUSSION 


Significance and Ranges of Parameters and 
Dimensionless Variables 


The subcooling parameter S = c p (T f - T c )/L 
appears because of the manner in which specific 
heat term was nondimens ionali zed in the basic 
equations. Although the numerator c-^Tf - T c ) is 
the maximum internal energy that coulcl be removed 
if a unit mass of the solidified material were sub- 
cooled from the freezing temperature to the lowest 
temperature in the system T c , the parameter does 
not give any indication of how much subcooling is 
actually experienced by the frozen layer. As the 
coolant temperature approaches absolute zero, S 
increases to a maximum value. For many materials 
employed in engineering practice such as molten 
metals and water, the maximum S is approxi- 
mately 3. 

The parameter R b (V k )/[ (1 / h c> + ( & Aw)] 
is a ratio of heat flow resistances. The Xg/k is 
the resistance offered by the steady-state thick- 
ness of the frozen layer, while (l/h c ) + (a/k^) is 
the resistance offered by the wall and the convec- 
tive coefficient on the coolant side of the wall. 
Since the thickness Xe can be very thin, the 
value of R can approach zero. For large Xg or 
if is large and a/ky small, R can be very 

large; hence, R can range from zero to infinity. 

To aid the clarity of the discussion and the 
physical understanding of the reader, it is con- 
venient to regard the parameter R as a measure 
of the steady-state frozen layer thickness Xg 
for fixed values of h c , a, k, and k w . 

The steady-state frozen layer thickness Xs 
in Equation (l) depends on the many independent 
variables in the system. For example increasing 
the liquid temperature Tj, the liquid heat trans- 
fer coefficient h}, or the coolant temperature 
T c ; or decreasing the coolant side convective coef- 
ficient h c all decrease the steady ice layer 
thickness. The dimensionless instantaneous solidi- 
fied thickness X' = X/X g varies by definition 
from zero to unity. 

The dimensionless time given by 

t M T Z - T f> 

pLXs 
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is the ratio of two heat flow quantities. The 
Th^ (Tj - T f ) is the total heat convected from the 
warm liquid to the surface of the frozen layer dur- 
ing the time t. The pLXs is the amount of 
latent heat that is removed to form the steady - 
state frozen layer. A large t* indicates that 
the heat convected to the frozen layer by the warm 
liquid has been large compared with the heat of 
solidification. 


Comparison of Solutions 

Figure 2 shows the growth of the solidified 
layer as predicted by four different equations pre- 
sented in the analysis for several values of the 
parameters R and S > 0. Also included are re- 
sults for the case where the subcooling in the 
frozen layer is neglected , S = 0 as computed from 
Equation (20). First consider the lower set of 
curves where R = 0. 01 which corresponds to thin 
frozen layers. In the limit when heat capacity is 
neglected, S = 0, all the solutions obtained from 
Equation (18) are in agreement, by definition, with 
that given by Equation (20). When S is increased 
to its maximum practical value of 3, all of the 
solution methods are in agreement and are only 
slightly removed from the S = 0 curve. 

The upper set of curves are for a thick frozen 
layer (R = 100). The curves for S = 3 represent 
the extreme case of large heat capacity effects 
coupled with a thick frozen layer and provide the 
maximum deviation between the calculation methods. 
One of the most significant findings is the very 
close agreement between the second and third 
analytical iterations (curves "AI.^" and "AI-^"). 

This indicates that the third iteration has con- 
verged very close to the final solution. This will 
also be shown by the temperature distributions 
which will be given in Figure 3. Hence for prac- 
tical purposes the third iteration can be regarded 
as the correct solution. The method using a pre- 
scribed temperature profile (PT curves) is in very 
good agreement with the analytical iterations, and 
hence the prescribed temperature method can also be 
used with good accuracy for the present problem. 

The interface method (IF curves) deviates somewhat 
from the other methods, but at a given x 1 still 
provides X' values within 5 percent of the cor- 
rect solution. The reason for this deviation 
stems from the fact that the growth times are re- 
lated to the first and second partial derivatives 
of the temperature at the interface (Eq. (41)). 

For the IF method these derivatives are obtained 
from the approximate temperature profile (Eq. (32)) 
rather than from a solution of the energy equation 
as in the AI method. 


Temperature Distributions Within Freezing Layer 

A few temperature distributions within the 
freezing layer are shown in Figure 3. Hie dimen- 
sionless temperature is shown as a function of po- 
sition within the layer for four different dimen- 
sionless layer thicknesses. An interesting char- 
acteristic is that all the temperature profiles are 
fairly close to being linear. The heat convected 
by the warm liquid to the interface combined with 
the latent heat of fusion that also arises at the 
interface, is conducted through the frozen layer 
and tends to establish a linear temperature pro- 
file. Hie heat capacity of the layer however tends 
to produce a curved profile. Hie fact that most of 


the temperature profiles are nearly linear suggests 
that the subcooling energy is small compared with 
the convective and latent heats at the frozen layer 
interface. 

Figure 3(a) shows the temperature distribu- 
tions as computed by various methods for a large 
value of both R and S. The rapid convergence of 
the successive analytical iteration solutions (AI) 
is demonstrated. The difference between the tem- 
perature distributions for the second and third 
iterations is very small. It appears that a fourth 
approximation would be indistinguishable from the 
third approximation and would have an insignificant 
effect on the frozen layer growth curves in Fig- 
ure 2. The profiles evaluated from the prescribed 
temperature (PT) and interface integration (IF) 
methods are not too far removed from the analytical 
iteration profiles. Figures 3(b) and (c) show how 
the agreement of all the methods is improved as 
either S or R is made smaller. 


Solidification Times as a Function of R and S 
Parameters 

While all the solutions agree well with the 
third analytical iteration which is taken to be the 
correct solution, they vary considerably in their 
form and difficulty of evaluation. The simplest 
solution is Equation (20) for the S = 0 condi- 
tion, but its range of application is restricted to 
small R values as indicated by Figure 2. The 
second analytical iteration. Equation (22), is 
quite accurate for all ranges of R and S and 
has a form that can be easily and quickly evaluated 
even on a desk calculator if necessary. The third 
AI, PT, and IF solutions all require a digital com- 
puter for their evaluation. For these reasons the 
second AI solution is the most convenient analyti- 
cal form given here for use in engineering applica- 
tions. The graphical results which will now be 
presented were evaluated from the third AI solu- 
tion (Eq. (26)). 

In Figure 4 is presented the dimensionless 
time t* against R with S as a parameter. 

Each group of curves gives the t 1 required to 
form a given frozen layer thickness X*. 

At each X’ it is interesting to note the be- 
havior of the group of curves as R becomes small. 
All the curves for S ^ 0 approach the curve for 
3=0, which means that the effect of heat capacity 
tends to become negligible as R diminishes. A 
small R means that the temperature drop across 
the frozen layer is small compared with the total 
drop across the wall resistance and convective re- 
sistance on the coolant side. This means that for 
small R the actual frozen layer subcooling is a 
small part of the maximum possible sub cooling 
Tf - T c . Hence the effect of the S parameter 
dies away at small R. 

As R becomes large it is seen that for any 
X', the t 1 is greater as S is increased from 0 
to 3. From this fact two important observations 
can be made concerning the subcooling of the frozen 
layer. The first observation is that when the heat 
capacity is not included in an analysis of a grow- 
ing frozen layer, serious errors can result in the 
predicted growth times. This conclusion will be 
shown as follows by comparing, for a fixed R, two 
cases, one where the specific heat is neglected so 
that S = 0, and another where S > 0. Under these 
conditions all the conditions such as the h^, h Q , 
Tj, T c , k, k^, Xg, remain the same for both cases. 
The ratios of the times t 1 and t for S = 0 
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and S > 0 are then equal ; that is. 



From Figure 4 it is found that the 
always greater than unity and hence 


x 1 ratio is 

V° > T Op=0' 


Thus, a simplified analysis neglecting heat capac- 
ity will always lead to predicting growth times x 
that are too short. 

The second observation is that, when the heat 
capacity is properly accounted for, the dimensional 
time t does not increase as x* increases when 
the S parameter of the frozen layer is increased. 
In fact the opposite is generally true, x de- 
creases. For a given solidifying material S is 
increased by lowering the coolant temperature T c . 

A lower coolant temperature, however, increases 
the solidification rate which more than compen- 
sates for the increased heat capacity that must be 
removed. This will now be discussed in more de- 
tail. 


Consider a fixed R and a fixed steady-state 
frozen layer thickness X s . For a given liquid 
where L, Cp, and Tf are assumed constant, an in- 
crease in S is accomplished by lowering the cool- 
ant temperature T c . Then since R depends on 
X g , and X g contains the ratio (T f - T c )/ 

hj(Tj - T f ), the quantity hj(Tj - T f ) must be pro- 
portionately increased to keep X s unchanged. If 
S is increased from to S 2 , then 


1 1 

i 

/-a 

E-^ 

r-> 

|2 

hj(Tj - T f )j 



Using the definition of the dimensionless time re- 
sults in 


t 2 

- T f 5 ’ 

j_2 

T 1 


|l 



1 


^2 

T 

1 


s 

s 


2 
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From Figure 4, it is found that for any X/Xg the 
ratio Sg/S-L is always larger than xg/x-^ so that 
ig < t^. The conclusion is that when a frozen 
layer of a specified steady- state thickness is be- 
ing formed in a given liquid, the layer will grow 
faster as the coolant temperature T c is decreased 
even though the subcooling of the frozen layer is 
thereby increased. 

Another feature of Figure 4 is that for any 
given X* the x' is essentially no longer de- 
pendent on R when R > 100. This implies that 
for large R's, the combined thermal resistance of 
the wall on which the layer is formed and the cool- 
ant convective resistance is negligible when com- 
pared to that of the frozen layer. As a conse- 
quence, the instantaneous frozen layer thickness 
is governed solely by the degree of subcooling im- 
posed on the layer. This is shown by the follow- 
ing simple relation derived from the second ap- 
proximate iterative solution (Eq. (22)) by taking 
a limit as R -* » 


3 + S 


II 


R-*» 


(-X- - ln|l - X'|) 


(43) 


where the expression in parentheses is the S = 0 
solution, Equation (20), for large R. A compari- 
son of Equation (43) with Figure 4 shows that the 
equation fits the curves quite well for R > 100. 

As a consequence, the factor (3 + S)/3 represents 
a correction factor that can be used to get rapid 
estimates from the S = 0 solution for R > 100 
and 0 < S < 3. 

One important aspect of Figure 4 is its use- 
fulness in making rapid estimates of the thickness 
X against time t for any given R and S. To 
illustrate, consider a given liquid, coolant, and 
wall material where the properties are known and 
the T^, T c , h^, h c are specified. With this 

information given, the quantities R and S can 
be calculated. Then from Figure 4, the x’ values 
can be found by simple interpolation for each X' 
shown. Having obtained values of t* against X', 
the relation between dimensional x and X can 
easily be calculated from the definitions of x' 
and X*. 


CONCLUDING REMARKS 


The goal of the analysis was to develop a 
means for predicting the transient growth of the 
frozen layer that forms when a flowing warm liquid 
is in contact with a cold flat plate that is con- 
vectively cooled on the opposite side. This goal 
was achieved by employing two basic approaches to 
the problem: (l) a growth relation obtained by in- 

tegrating the transient conduction equation over 
the entire frozen layer thickness and (2) a heat 
balance at the frozen layer- liquid interface. 

By use of a proper dimensionless layer thick- 
ness and time variable, the frozen layer growth 
can be expressed as a function of two parameters. 

One parameter provides a measure of the maximum 
possible subcooling energy of the frozen layer as 
compared with the latent heat of fusion. The sec- 
ond is the ratio of the heat flow resistance of the 
steady-state layer, to the combined resistance of 
the wall and convection coefficient on the coolant 
side. 

From the two analytical approaches five solu- 
tions were developed, four from the first approach 
and one from the second. All solutions except one 
gave accurate predictions of thickness against time 
for all values of the parameters. The one excep- 
tion was valid for only thin frozen layers where 
subcooling was unimportant, as this solution had 
neglected the heat capacity of the frozen layer. 

Of the accurate solutions, two that were analyti- 
cally derived successive iterations were the most 
convenient to use and the most precise. The form 
of one of these iterative approximations was suf- 
ficiently simple that it could easily be evaluated 
on a desk calculator. The other, from which the 
graphical results were prepared, appeared well con- 
verged to the correct solution. 

The results of this analysis led to several 
conclusions. The first is that the sub cooling of a 
frozen layer during its formation slows the growth 
rate substantially for some conditions. This is 
particularly true of thick layers. The second is 
that the temperature profiles in the frozen layer 
were found to be almost linear at all times. This 
suggests that the energy required to subcool the 
frozen layer is much less than the combined convec- 
tive and latent heat passing through the layer. 
Lastly, the analytical iterative technique developed 
here resulted in a rapidly converging means for 
solving the nonlinear transient freezing problem. 
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NOMENCLATURE 

A coefficient of linear term in Eq. (32) 

a thickness of cooled vail 

B coefficient of quadratic term in Eq. (32) 

c p specific heat of solidified material 

h convective heat transfer coefficient 

k thermal conductivity of solidified material 

k w thermal conductivity of vail 

L latent heat of fusion 

q heat flux in x direction 

R dimensionless parameter: 

(X B /k)/[(l/h c ) + (a/kj] 

S dimensionless parameter: c p (Tf - T c )/L 

T temperature 

T’ dimensionless temperature, 

(T - T f )/(T C - T f ) 

X thickness of frozen layer 

X' dimensionless thickness of frozen 

layer, X/X s 

X g thickness of frozen layer at steady state 

x position coordinate in frozen layer 

x' dimensionless coordinate, x/Xs 

t dimensionless coordinate, x/X 

p density of solidified material 

x time 

x* dimensionless time, xh^Tj - T f )/pLXg 

Subscripts : 

c refers to coolant 

f at freezing point 

1 liquid phase of solidifying substance 

2 interface betveen frozen layer and vail 
I, II, successive iterative approximations 

III 
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Flowing coolant 
with constant 



Figure 1. - One-dimensional model of transient solidification of 
a warm liquid on a convectiveiy cooied piaie. 
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Dimensionless time, t' = xh^Tf- Tf)/pLX s 

Figure 2. - Comparison of methods for predicting the instantaneous thickness 
of the frozen layer. (Some curves dotted for clarity only. ) 









Dimensionless time, 



(a) S = 0, 0.5, and 1. 


gure 4. - Dimensionless growth times against the thermal resistance parameter for various dimensionless thicknesses and subcooling 
parameters. 




